function [xval, f1, f2, err_a, err_b, b_err_a, b_err_b, post_1_a, post_1_b, post_2_a, post_2_b evi_a, evi_b, cond_err_a, cond_err_b, boundary_err_a, boundary_err_b, bound_a, bound_b] = bsp3_2()
%% Die Funktion bsp3_2.m erstellt 2 Klassen mit zufaelligen normalverteilten Featurevektoren. 
% Es werden fuer 2 verschiedene a priori Wahrscheinlichkeiten (P(w1)=0,5
% und P(w2)=0,5 bzw. P(w1)=0,9 und P(w2)=0,1)
% die Randverteilungen und Posteriors, sowie im folgenden angefuehrten Fehlerwahrscheinlichkeiten und
% Entscheidungsgrenzen berechnet
% xval - Werte von -3 bis +9 in 0.1 Schritten
% f1, f2 - Werte der Normalverteilungen
% err_a, err_b - Fehlerraten
% b_err_a, b_err_b - Bayes Fehlerraten
% post_1_a, post_2_a, post_1_b, post_2_b - Posteriori Wahrscheinlichkeiten
% evi_a, evi_b - Randwahrscheinlichkeiten
% cond_err_a, cond_err_b - bedingte Fehlerwahrscheinlichkeiten
% boundary_err_a, boundary_err_b - Fehlerraten bei Entscheidungsgrenze xp=4
% bound_a, bound_b - Bayes-Entscheidungsgrenzen
% das suffix _a ist bezeichnend fuer die Aufgabenstellung 
% P(w1)=0,5 und P(w2)=0,5
% das suffix _b ist bezeichnend fuer die Aufgabenstellung 
% P(w1)=0,9 und P(w2)=0,1

%% 2 Klassen mit zufaelligen normalverteilten Featurevektoren erstellen
x=[-3:0.1:9];
x1 = normpdf(x, 1.5, 1); % Normalverteilung mit mu = 1.5 und sigma = 1
x2 = normpdf(x, 3, 1); % Normalverteilung mit mu = 3 und sigma = 1

% Featurewerte aufsummieren
ts1 = sum(x1);
ts2 = sum(x2);

% Wahrscheinlichkeiten berechnen
pxw1 = x1 / ts1;
pxw2 = x2 / ts2;

%% Randverteilung fuer P(w1) = P(w2) = 0.5
pc1 = 0.5;  % Wkt fuer Klasse 1 = 0.5
pc2 = 0.5;  % Wkt fuer Klasse 2 = 0.5

pxw1a = pxw1 * pc1;
pxw2a = pxw2 * pc2;

pxwra = pxw1a + pxw2a;

% disp('Tabelle: Randverteilung - bei P(w1) = P(w2) = 0.5');
randverteilung = [pxw1a, sum(pxw1a) ; pxw2a, sum(pxw2a); pxwra, sum(pxwra)];

%% A posteriori Wahrscheinlichkeiten
temp_a = pc1 * pxw1 + pc2 * pxw2;

for i = 1 : size(x, 2);
post1a(i) = (pc1 * pxw1(i)) / temp_a(i);
post2a(i) = (pc2 * pxw2(i)) / temp_a(i);
end

%Zeilenvektoren zu Spaltenvektoren transformieren
post1a = post1a';
post2a = post2a';

% A Posteriori Wahrscheinlichkeiten ausgeben
% disp('Tabelle: A Posteriori Wahrscheinlichkeiten - bei P(w1) = P(w2) = 0.5');
a_posteriori_wkt = [post1a, post2a];

% Fehlerrate berechnet sich aus der Summe der zeilenweisen multiplikation
% der geringeren a posteriori Wkt der beiden Klassen mit Nenner aus dem
% Bayes Theorem
error_rate_a = 0;
for i = 1 : size(x,2)
    tempsum = temp_a(i) * min(post1a(i), post2a(i));
    error_rate_a = error_rate_a + tempsum;
end

disp('Error-Rate - bei P(w1) = P(w2) = 0.5');
error_rate_a;


%% Randverteilung, a Posteriori WKT mit P(w1) != P(w2)

% Randverteilung fuer P(w1) = 0.9, P(w2) = 0.1
pc1 = 0.9;  % Wkt fuer Klasse 1 = 0.9
pc2 = 0.1;  % Wkt fuer Klasse 2 = 0.1

pxw1b = pxw1 * pc1;
pxw2b = pxw2 * pc2;

pxwrb = pxw1b + pxw2b;

% disp('Tabelle: Randverteilung - bei P(w1) = 0.9, P(w2) = 0.1');
randverteilung = [pxw1b, sum(pxw1b) ; pxw2b, sum(pxw2b); pxwrb, sum(pxwrb)];

% A posteriori Wahrscheinlichkeiten
temp_b = pc1 * pxw1 + pc2 * pxw2;

for i = 1 : size(x, 2);
post1b(i) = (pc1 * pxw1(i)) / temp_b(i);
post2b(i) = (pc2 * pxw2(i)) / temp_b(i);
end

%Zeilenvektoren zu Spaltenvektoren transformieren
post1b = post1b';
post2b = post2b';

% A Posteriori Wahrscheinlichkeiten ausgeben
%disp('Tabelle: A Posteriori Wahrscheinlichkeiten - bei P(w1) = 0.9, P(w2) = 0.1');
a_posteriori_wkt = [post1b, post2b];

% Fehlerrate berechnet sich aus der Summe der zeilenweisen multiplikation
% der geringeren a posteriori Wkt der beiden Klassen mit Nenner aus dem
% Bayes Theorem
error_rate_b = 0;
for i = 1 : size(x,2)
    tempsum = temp_b(i) * min(post1b(i), post2b(i));
    error_rate_b = error_rate_b + tempsum;
end

% disp('Error-Rate - bei P(w1) = 0.9, P(w2) = 0.1');
error_rate_b;

%% Bayes Conditional error
bayes_cond_err_a = min(post1a, post2a);
bayes_cond_err_b = min(post1b, post2b);

%% Bayes Error rate - 
% errechnet sich aus dem Integral ueber x von -inf bis
% +inf der Funktion (bedingter Bayes-Fehler) * (Randverteilung)
% in unserem Fall reicht der Bereich von -10 bis +20 (bzw von -3 bis 9,
% jedoch aus Testgruenden auf -10 - +20 belassen - nicht 100% funktionstuechtig)
y = quad(@integrand_a, -10, 20);

  function value = integrand_a(x)
      for i=1:size(x,2)
    value(i) = bayes_cond_err_a(i)*pxwra(i); 
      end
  end

bayes_error_rate_a = y;

y = quad(@integrand_b, -10, 20);

  function value = integrand_b(x)
  for i=1:size(x,2)
    value(i) = bayes_cond_err_b(i)*pxwrb(i); 
      end   
  end

bayes_error_rate_b = y;

%% Conditional Error bei der Grenze x=4
% Entspricht den posteriori Wktn der Klasse 2 fuer x < 4, ansonsten der
% posteriori Wktn der Klasse 1
index = 0;
for t=-3:0.1:9
    index = index + 1;
if( 4 > t)
    cond_error_a(index) = post2a(index);
else
    cond_error_a(index) = post1a(index);
end
end

index = 0;
for t=-3:0.1:9
    index = index + 1;
if( 4 > t)
    cond_error_b(index) = post2b(index);
else
    cond_error_b(index) = post1b(index);
end
end

%% Error Rate bei Grenze x=4 
% errechnet sich aus dem Integral ueber x von -inf bis
% +inf der Funktion (bedingter Bayes-Fehler) * (Randverteilung)
% in unserem Fall reicht der Bereich von -3 bis 9 - nicht 100%
% funktionstuechtig
y = quad(@integrand_c, -3, 9);

  function value = integrand_c(x)
  for i=1:size(x,2)
    value(i) = cond_error_a(i)*pxwra(i); 
  end
  end

b_error_rate_a = y;

y = quad(@integrand_d, -3, 9);

  function value = integrand_d(x)
  for i=1:size(x,2)
    value(i) = cond_error_b(i)*pxwrb(i); 
  end
  end

b_error_rate_b = y;


%% "Grafisches" ermitteln der Entscheidungsgrenze
% Betrachtet man die Graphen der Posteriori Wktn so befindet sich die
% Entscheidungsgrenze bei dem x-wert, an dem sich der Schnittpunkt zwischen
% den beiden Kurven befindet
for i=1:size(x,2)-1
    if(post1a(i) < post2a(i))
        boundary_a = x(i)
        break;
    end
end

for i=1:size(x,2)-1
    if(post1b(i) < post2b(i))
        boundary_b = x(i)
        break;
    end
end

    
%% Rueckgabewerte setzen
xval = x;
err_a = error_rate_a;
err_b = error_rate_b;
b_err_a = bayes_error_rate_a;
b_err_b = bayes_error_rate_b;
post_1_a = post1a;
post_1_b = post1b;
post_2_a = post2a;
post_2_b = post2b;
evi_a = pxwra;
evi_b = pxwrb;
cond_err_a = cond_error_a;
cond_err_b = cond_error_b;
boundary_err_a = b_error_rate_a;
boundary_err_b = b_error_rate_b;
f1 = x1;
f2 = x2;
bound_a = boundary_a;
bound_b = boundary_b;

    
end
